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The magnetic properties of iron (spin and orbital magnetic moments, magnetocrystalline 
anisotropy energy) in various geometries and dimensionalities are investigated by using a 
parametrized tight-binding model in an s, p and d atomic orbital basis set including spin polar- 
ization and the effect of spin-orbit coupling. The validity of this model is well established by 
comparing the results with those obtained by using an ab-initio code. This model is applied to the 
study of iron in bulk bcc and fee phases, (110) and (001) surfaces and to the monatomic wire, at 
several interatomic distances. New results are derived. The variation of the component of the orbital 
magnetic moment on the spin quantization axis has been studied as a function of depth, revealing a 
significant enhancement in the first two layers, especially for the (001) surface. It is found that the 
magnetic anisotropy energy is drastically increased in the wire and can reach several meV. This is 
also true for the orbital moment, which in addition is highly anisotropic. Furthermore it is shown 
that when the spin quantization axis is neither parallel nor perpendicular to the wire the average 
orbital moment is not aligned with the spin quantization axis. At equilibrium distance the easy 
magnetization axis is along the wire but switches to the perpendicular direction under compression. 
The success of this model opens up the possibility of obtaining accurate results on other elements 
and systems with much more complex geometries. 

PACS numbers: 71.15.Ap,71.20.Bc,71.70.Ej,73.20.At,73.22.Dj,75.10.Lp,75.70.Rf,75.50.Bb,75.75.+a 

I. INTRODUCTION 

The magnetic properties of nanoparticlcs. thin films and wires have attracted recently a lot of attention due to their 
potential technological applications. It is thus very important to investigate the influence of dimensionality on these 
properties. In many experimental systems, some atoms are in a bulk- like environment while some others have a very 
low coordinence in a strongly asymmetric environment. This is the case of clusters, at the edge of adsorbed islands, 
supported wires, along step edges or in the nanoconstriction region of a break junction etc.. As a consequence the 
magnetic properties of such systems need to be treated at an atomic level, by studying their electronic structure in 
the framework of quantum mechanics. 

In principle these properties can be determined from ab-initio calculations. However the computer time and storage 
increase drastically with the number of inequivalcnt atoms. Thus there is a need for simplified methods still based on 
quantum mechanics which capture the essential physics. In this context tight-binding (TB) methods are ideally suited 
for these calculations. Indeed TB methods can be extended to magnetic systems by introducing Hubbard two-body 
terms treated in the Hartree-Fock approximation , and can easily include spin-orbit effects 0|. Moreover they 
allow the calculation of local physical quantities (spin or orbital moment etc..) in a straightforward manner which 
gives a physically transparent understanding of the phenomena. 



2 



The anisotropy of the magneto-crystalhne energy (MAE), although small in transition metals, is of fundamental 
importance since it determines the easy magnetization axis. The MAE results from the coupling of the spin and 
orbital moments. In the bulk it is well known that the orbital moment is nearly quenched due to the high symmetry 
of the potential. When dimensionality or symmetry is reduced the orbital moment is less and less quenched and 
the MAE increases rapidly. The MAE is routinely measured by magnetic hysteresis or torque measurements 
1^, for instance. More recently the development of X-ray magnetic circular dichroism techniques has allowed the 
experimental determination of the orbital moment Q- On the theoretical side the calculation of the MAE in bulk 
pure ferromagnetic transition metals is a challenge because of its minuteness (typically some ^eV). The first attempts 
to determine the MAE on surfaces have been performed on free standing mono-layers either using a perturbative 
treatment of the spin-orbit coupling in a tight-binding model^, or with self-consistent ab-initio technique^, 
[lo| |. The calculations performed by Bruno [a, Ql on mono-layers have been extended to study slabs with several 
atomic layers^^ in a pure d-band model. More recently the case of layered ordered alloys of the CuAu or 
CsCl types containing at least one ferromagnetic element ^ | , deposited ferromagnetic over-layers on non-magnetic 
substrates [l^. and multilayers have been investigated [15I 16 1. Simultaneously the orbital magnetism driven by the 

H'n-orbit interaction has been computed using the tight-binding model (perturbatively0,y| or non Derturbativelvjirl 
) as well as ab-initio codes[l9|,|20j. 

The aim of this paper is to develop a TB model allowing the determination of the magnetic properties of transition 
metals in various geometrical configurations, going from highly coordinated and symmetric environments to low 
coordinated and anisotropic geometries. In a first step we have found useful to check the validity of the model by 
a detailed comparison with the results provided by ab-initio methods on simple systems. In addition, the use of 
both methods yields a better physical understanding of these properties, and, moreover, the simplicity of our model 
allows a more detailed analysis of each system. We use a non-orthogonal basis set of s, p and d valence orbitals. The 
parametrization of the non-magnetic and non-relativistic Hamiltonian was derived by Mehl and Papaconstantopoulos 
[2l|. The possibility of spin polarization is then introduced using a Stoner like model in which the splitting between 
the energy levels of up and down spin orbitals is governed by the Stoner parameter Im' and is proportional to the 
magnetic moment carried by d electrons. Indeed it is well known that s and p electrons are very weakly polarized. 
The relativistic effects, limited to spin-orbit coupling between d electrons, are taken into account by adding the intra 
atomic matrix elements of this coupling determined by a parameter ^. The two parameters Idd' and ^ are fixed by 
comparison with ab-initio calculations. For the Stoner parameter the variation of the magnetization as a function of 
interatomic distance in the bulk phase can be used. The determination of the spin-orbit coupling parameter relies on 
the study of the degeneracy removal of energy bands at high symmetry points or directions of the Brillouin zone. 

We have applied our model to the study of iron in various atomic arrangements going from the bulk, in the exper- 
imentally observed phases bee and fee, to simple surfaces and finally to the monoatomic wire, at several interatomic 
distances. It is found that with a unique value of the Stoner parameter we are able to reproduce the variation of the 
spin magnetic moment in a wide range of lattice spacings in the bee bulk phase. It is well known that the magnetic 
properties of bulk fee iron arc rather complex. In spite of this complexity the calculations carried out in this work 
are in excellent agreement with ab-initio predictions. The spin-orbit coupling parameter is then determined, as ex- 
plained above, and a single value of this parameter is able to reproduce the ab-initio band structure. Furthermore the 
contribution of the orbital moment to the magnetization in bee Fe is very close to the experimental value. The case 
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of surfaces represents a more stringent check since the dimensionahty is reduced and aU atoms are not geometrically 
equivalent. It is in particular interesting to follow the variation of the magnetic properties when going from the 
outermost to inner (bulk-like) layers. The spin-polarized surface projected band structure as well as the variation 
of the spin magnetic moment as a function of depth derived from ab-initio calculations are perfectly reproduced by 
our model. In particular the (001) surface atoms have a saturated moment contrary to those of the (110) surface. 
Introducing the spin-orbit coupling in the TB Hamiltonian allows the calculation of the MAE and of the orbital 
magnetic moment at each inequivalent site. Due to the efficiency of our model it is possible to check rigorously the 
convergence of these quantities when increasing the number of k points. The orbital moment is strongly increased 
on surface atoms (about twice the bulk value for the (001) surface), and recovers its bulk value on the third and 
innermost layers. The wire is the most anisotropic atomic arrangement with the lowest coordinence. Even though 
the TB parameters are fitted on bulk ab-initio data only, the non-magnetic band structure calculated with these 
parameters is in good agreement with ab-initio calculations, proving their good transferability. The two additional 
parameters I^d' and ^ are also perfectly transferable. Indeed, without changing the Stoner parameter, the variation 
of the spin magnetic moment with the interatomic distance is satisfactorily reproduced (in particular the saturated 
solution appears abruptly at the same interatomic spacing) and the splitting of bands due to spin-orbit coupling is 
exactly the same, compared to ab-initio results. The calculation of the MAE reveals that at theoretical equilibrium 
the easy axis is parallel to the wire, but at smaller interatomic distances corresponding to unsaturated magnetic 
solutions the easy axis is perpendicular to the wire. We have finally checked the validity of Bruno formula 0| relating 
the MAE to the anisotropy of the orbital moment, and found that this relation is almost strictly obeyed around the 
equilibrium distance. Indeed at this interatomic distance the up spin bands are filled and the exchange splitting is 
large compared to the d bandwidth, which was not the case for the (001) surface that, although saturated, has much 
a wider d bandwidth. 

The paper is organized as follows. In section 2 the formalism of our model is presented in details, in particular 
the derivation of the x, y and z components of the orbital and spin moment formula in a non-orthogonal basis set. 
The spin-orbit coupling being small we have recalled the perturbation treatment of the MAE and of the orbital 
moment, from which the analytical expression of the anisotropy laws are directly obtained, and arc used to analyze 
our numerical results. In section 3 the Stoner parameter is determined and used to study in detail the magnetic 
properties of bee and fee iron. Section 4 is devoted to the study of (001) and (110) surfaces. Finally in section 5 we 
present an exhaustive study of the monoatomic wire. Conclusions are drawn in section 6. 

II. FORMALISM 

A. Spin polarized tight-binding model 

We choose as a basis set the real s, p and d valence atomic orbitals centered on each site i. They are denoted 
by A and fi indices (A,/i — 1,9) and numbered as follows: s,pxTPy,Pz,dxy,dyz,dzx,dx2-y2,djz^_r2, the x,y,z coor- 
dinates being taken along the crystal axes. The tight-binding (TB) hamiltonian for the non-magnetic (NM) state is 
then completely determined by its intra-atomic matrix elements (i.e., the s, p and d atomic levels) ex and its inter- 
atomic matrix elements (i.e., the hopping integrals) f3^l^ which have been tabulated as a function of 10 Slater-Koster 
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(SK) parameters (ssa, spa, sda, ppa, ppir , pda, pdi: , dda, ddir, ddS) and of the direction cosines of the bonding direction 
Ryl^l- Following the (MP) scheme developed by Mehl and Papaconstantopoulos [2^, the atomic levels depend on 
the atomic environment (number of neighbors and interatomic distances) while the SK parameters are function of Rij 
only. Finally the Schroedinger equation in the atomic orbital basis involves also overlap integrals S^J^ depending on 
the bonding direction Rij when the non-orthogonality of the basis set is taken into account. All these quantities {ex, 
Pij^, S^j^) are written as analytic functions depending on a number of parameters which are determined by a least 
mean square fit of the results of ab-initio electronic structure (band structure and total energy) calculations either 
in the Local Density (LDA) or in the Generalized Gradient (GGA) approximations. These parametrizations will be 
denoted as TBLDA and TBGGA in the following. The analytical form of the functions can be found in Ref.|2^ and 
the numerical values of the parameters for Fe can be found in Ref. |23| . 

n 

In order to account for spin polarization, we use a simplified Hartree-Fock (HF) scheme [l] to define atomic levels 
depending on spin (eao-)- These diagonal elements of the hamiltonian can be written in the basis of spin-orbitals jiAcr) 
in which the spin quantization axis is parallel to the magnetization (a = -f 1(— 1) for up (down) spin). When all atoms 
in the system are geometrically equivalent we get: 

£s,a = £os + Uss^ + {Usj, - ^)Np + [Usd - ^)Nd 
(UssMs + JspMp + JsdMd) 
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£pQ(T = £op + UspNg + {Upp' 2~'^^^ ~ 2^^^P' ~ '^'^PP'^^Pc. + UpdNd 

^ 2"(^sp^^s + Jpp'Mp + JpdMd + {Upp' + Jpp')mp^) 
£d^a = £od + {Usd - ^)Ns + {Upd - ^)Np + {Udd' - ^)Nd - ^{Udd' - SJdd')nd^ 

- ^{JsdMs + JpdMp + Jdd'Md + {Udd' + Jdd' )md^ ) (1) 

where Ng(^p ,i-^ and Ms(^p^d) ar^ the total number of electrons and the total moment on each atom, respectively, in s, 
p and d orbitals while 'T.p^(d„) a-nd ™p„(dQ) denote the total occupation number (i.e., for both spins) and moment 
in orbital Pa{da). Finally the U and J parameters arc Coulomb and exchange integrals which involve two different 
orbitals, save for Uss, and are assumed to depend only on the orbital quantum numbers of these orbitals. 

We further assume that the asphericity of both the charge distribution and magnetic polarization can be neglected, 
i.e., Up^ = Np/S, Ud^ ~ Nd/5i, mp^ = Mp/3, rud^ = Md/5. In these conditions, when the system is non magnetic 
all the non vanishing terms in Eqs^arc accounted for implicitly by the expression of £\ in the MP scheme |2ll |. In 
addition the equations giving the atomic levels can be further simplified by noting that the spin polarization of s and 
p electrons is very small. As a consequence Eq^can be approximated by: 





= £s 


- JsdMd 






a 




= £p 


— JpdMd 










= £d 


- -^Idd'Md 



where es{p.d) are the NM levels j2l| and Idd' = {Udd' + (>Jdd')/5 can be identified with the Stoner parameter. The 
numerical value of this parameter will be determined in section 3.1 in order to reproduce as closely as possible the 
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variation of the magnetic moment as a function of the interatomic distance that can be obtained from an ab-initio 
calculation. Finally, Jsd and Jpd arc one order of magnitude smaller than Idd' and we have taken Jsd — Jpd ~ 
Idd'/^O- This completely defines our TB spin polarized hamiltonian Htbhf hi the absence of spin-orbit coupling for 
a system of equivalent atoms, and when the overlaps are neglected. 

In the general case where overlaps are taken into account and all atoms in the systems are not geometrically 
equivalent, the Hamiltonian becomes: 



H'^i-^^ are the matrix elements provided by the MP parametrization of the Hamiltonian. The second term in which 
5Ni is the net total charge on atom i, prevents large charge transfers when inequivalent atoms are present, and will be 
discussed in section 4. Finally in the last term which accounts for spin polarization A^^ = JsdM^, JpdAI^ and Idd'M^ 
for s,p,d orbitals, respectively. 



B. The spin-orbit coupling 



The spin-orbit interaction for a single atom is given by: 



Hso^-^{VVAp).a- (4) 



where V is the atomic potential, p is the momentum operator and cr are the Pauli matrices. Taking into account the 
spherical symmetry of the potential, Hgo can be rewritten as: 



(5) 



with: 



2w?(? r dr 

L = r A p and S = her /2 are, respectively, the angular orbital and spin momentum operators. The matrix elements 
of Hso in the basis of atomic spin-orbitals |Acr) are: 



(Aa|H,o|M^') =eAM(AfT|L.S|/la' 



(7) 



with: 



Cam^ / 7^A(r)7^^(r)^(r)r^dr 

"'0 



(8) 
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where TZ\{r) is the radial part of the atomic orbital A and A denotes its angular part. Since ^(r) is well localized 
around r = 0, ^Xfi has a non negligible value only when TZ\{r) and 7?,^(r) are also well localized, i.e., for transition 
metals, when both A and fi are d orbitals in which case TZ\{r) = TZ^ir) ~ TZdir) and ^\f^ = ^ > 0). 

In the tight-binding approximation the crystal potential is written as V(r) = ^(k ^ I^il) ^^'^ ^so becomes: 



where is the angular orbital momentum operator with respect to the center i. For transition metals and due to the 
localized character of ^(|r — Ri|) we can neglect all matrix elements of Hgo save for the intra-atomic ones between d 
orbitals. These matrix elements are the same at each site and given in Appendix A in the spin framework a;", y", 2". 
In this framework 2" is the spin quantization axis defined by its polar and azimuthal angles 9, ip relative to the crystal 
axes. The x" and y" axes have been chosen in the following way: the x, y axes of the crystal are first rotated by the 
angle (p around z, this gives a new framework x', z' which is then rotated by an angle 9 around y' . The orbital and 
spin moments are usually expressed in units of h so that ^ is a parameter which has the dimension of an energy. Its 
numerical value will be deduced from ab-initio calculations in the following. 

C. Determination of the components of the spin and orbital moments in the spin framework. 



Let us first compute the average value of the three components of the total spin < Sx" > , < Syn > , < Sz" > in the 
spin framework. If we choose as a basis set of spin-orbitals the direct product of the orbitals |iA) with the eigenvectors 
of the operator Sz" denoted as t and j, the electron eigenfunctions lipn) in the crystal can be written: 



(Note that in the absence of spin-orbit coupling, there is no spin mixing in these eigenstates and, since the matrix 
elements of Htbhf are real, it is always possible to find a set of eigenvectors whose components are real and denoted 
as cf^^ in the following). The average values of the three spin components are given by: 



//so = ^e(|r-R.|)Li.S 



1. Spin moment. 




< s > 



n occ 



in a non orthogonal orbital basis set, we obtain: 



< Sx" > = 



n occ 




71 OCC 




(9) 



n occ 
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i.e., in the absence of spin-orbit coupling < Sx" >=< Sy" >= 0. 

For a system with full translational symmetry and a single atom per unit cell, the Bloch theorem yields (cr =t or 
i): 



^ :exp(ik.ROc$:,(k) (10) 



since each eigenstate n is labelled by a band index (a = 1, 9) and a wave vector k. iVat is the number of atoms. The 
spin components are the same on all sites i and arc given by: 



<Sx„>^Re J2 cXT*(kKiW5AM(k) 

A/j. 

(Q,k) occ 

<V>-Im E c^T*Wi(k)5AM(k) (11) 

Am 
(a,k) occ 

Xfia 
(a,k) occ 

with 5A^(k) = TVa^i E,, cxp(*k.(R, - R.))4^ 

When all atoms are not geometrically equivalent, we can define a spin on site i by identifying in EqsElall the terms 
involving this site and, similarly to what is done to define Mulliken charges, the overlap cross terms (i.e., those in 
which only one of the site indices is equal to i) are multiplied by a factor 1/2 to avoid a double counting of these 
terms. This condition ensures that these "local" spins arc real and that their sum is equal to the total spin. For 
example in a periodic system with several atoms per unit cell the local spin on each atom in the cell are given by 
equations similar to equation II II with an additional index, labelling the atom in the cell. For instance in a slab the 
local spin moment < Saz" > on layer a is given by: 

<^-">=i( E '^(cS.(k!t)c?M-(k|t)5^;:(k|i)+c?;,(ki|)c:;,,(k||)5^"^(k||))] (12) 

\ bXfia / 
(a,k|| ) occ 

with 5^1^ (k||) = ^||at Eij "^^pI*^!! -(-f^J ^ '^i))'^ia]b- ^l|at is the number of atoms in each layer of the slab and k|| the 
wave vector parallel to the surface, each atom being now labelled by a cell index, i or j, and a layer index, a or b. 
Corresponding changes must be made for the two other components of S. Finally let us recall that the spin magnetic 
moment M is related to the spin Sby<M>=— 2<S> (in Bohr magnetons iib)- 

2. Orbital moment. 

Up to now the orbital moment in the TB approximation has always been calculated by assuming an orthogonal 
basis set of atomic orbitals and only its z" component was determined. In these conditions, the component of the 



local orbital moment on site i in this direction is usually written as 



0: 



<U,„>^^mj p,Una{E)dE (13) 

where piima{E) is the local density of states at site i projected on the atomic orbitals \ilm) = TZi{r")Yi„i{6" , if") and 
spin function a, the variable r", 9" , Lp" being spherical coordinates relative to the spin framework, i.e.: 

Pama{E) = ^{iP„\ilma){ilma\tJjn)S{E - E„). (14) 

hn 
n 

Thus 

< Liz'i >= {■4;n\ilm(j)m{ilm(j\tpn) (15) 

Irna 
n occ 

This defines the operator Liz" which is diagonal in the \ilma) basis. Eg 1 151 can be generalized for the two other 
components of the orbital moment by noting that the corresponding operators are not diagonal in this basis. This 
gives: 

<L- >= ^ {ipn\ilma)[L^]i„i,i',n'{'il'm'a\'iljn) (16) 

n occ 

with L- = [Lix" , Liyff ^Liz"). Finally in the basis of real orbitals |zA(t) defined in the crystal frame, we have: 

<L;' 5] c^[L:']vc^,, (17) 

n occ 

The operators can be expressed as a function of the three operators Li^^ Liy, Liz projecting the orbital moment 
on the crystal axes, i.e.: 

Lix" = cos 6 cos ip Lix + cos 9 sin p Liy — sin 9 Liz 
Liyii = — sinp Lix + cos p Liy 

Liz" ~ &\n9 cosp Lix + wiiO smifi Liy + cos9 Liz (18) 

and the matrix elements of between two atomic orbitals A and p, centered on atom i defined with respect to the 
crystal axes are easily calculated (see Appendix A). These matrix elements are either vanishing or imaginary, thus 
[Li]A,j = — [Li]^A- In the absence of spin-orbit coupling, as stated above, the coefficients c^i"^^ are real and the orbital 
moment vanishes, i.e., < >= 0. 

EqElcan be generalized to take overlap into account (see Appendix B). This yields 

< L;' >= Re ^ c^:^[L:'],,5^''c^, (19) 

n occ 
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and. for a system with a full translational symmetry and a single atom per unit cell, this gives using Ea llOl 

<L:'>=Re J2 c^;(k)[L:'];,^5^.(k)c^,(k) (20) 

ctk occ 

This latter equation can be easily generalized to periodic systems with several atoms per unit cell, in the same way 
as for the spin moment. Finally, let us note that in the presence of spin-orbit coupling, the direction of the total 
magnetization (— < L + 2S >) may not strictly be parallel to the spin quantization axis z" . However Hgo being a small 
perturbation in Fe, we will often denote the spin quantization axis as the magnetization direction in the following. 

D. Magnetocrystalline anisotropy and orbital moment from perturbation theory 

We have seen in section 2.2 that spin-orbit effects can be limited to d orbitals. Furthermore the overlaps between 
these orbitals are close to zero and the spin-orbit coupling is a weak perturbation since ^ is much smaller than the Fe 
d bandwidth. Consequently spin-orbit coupling effects can be understood using a simple perturbation theory with a 
basis set of orthogonal d orbitals 

Let us consider the perturbation of the total energy due to Hso- Since the matrix elements of Hgo are a function 
of 9 and ip, this introduces an angular dependence of this perturbation which is known as the magnetocrystalline 
anisotropy. The first order term can be written: 

{n(j\Hso\na) (21) 

where \na) is an unperturbed state of energy E^^, i.e., 



ncr occ 



K>=E^"A.NAa) (22) 

Thus: 



AE^'^ = Cj2C>^a\L.Sma) ^ c°r.C (23) 

Xfi i 

na occ 

It is easily seen that lS.E^^^ vanishes since for each spin the (5 x 5) matrix {Xa\L.S\flcF) is imaginary (see appendix A). 
The second order perturbation of the total energy is given by: 



na occ 
n' a' unocc 



This yields: 



AE^'^^-e E T.C^<T\L.S\fla'){^I'a'\L.S\X'a)Y,h,iX,X\^,',^L,a,a') (25) 

XfiX' fj,' aa' ij 
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with: 



and: 



pO/,^' (E) = ^ c^LcTx'JiE - Si) (27) 



By using the relations between the matrix elements of L.S shown by Bruno [6|,|3| and recalled in Appendix A, Ai?*^^^ 
can be rewritten as: 

Ai?^^-* = isotropic term 

- e J2 (AT |L.S|/iT)(M'T |L.S|A't) '^'^%(A,A',/i',//,a,a') 

XfiX' fj,' ij.aa' 

= 5]Ai?f) (28) 

i 

where /S.Ef^ is the contribution of atom i to the perturbation energy. 

In the case of a system with full translational symmetry and a single atom per unit cell and using Ea llOl Eg 1281 can 
be transformed into 



A£;(2) = isotropic term - $2 ^ (A j |L.S|^ j) j |L.S|A' t) 



X 

k 



— oo 



with: 



Mxy{KE)^Mxx']iKE)-Mxx'i{KE) (30) 

and: 



Afxx'AKE) = Y,clZ*ik)cl?AmE - E^^A^)) (31) 

a 

the superscript refers to the unperturbed state as above and i?Kk) are the unperturbed eigenenergies. These 
equations can be generalized to a periodic system with several atoms per unit cell |ll| . 

Let us now consider the projection of the orbital moment on the spin framework axes. From Ea llTl it can be seen 
that the operators associated with these projections at a given site i can be written in an orthogonal basis set: 



L; =Y,\^Xa)[L;]x^{i^^a\ (32) 

Xfia 
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Within perturbation theory, we have: 



{n(T\L^\n'cT'){n'(T'\H^a\n(T) 

>--- m — TT^o + ^-^-^i 



na occ riV 



n (T unocc 



By substituting Easl32l and 1221 for and \na)^ respectively, into the preceding equation, we get: 

<L:'>=-2e E<^^|Ll'lMfT>(-""'^|L-S|A'a)^/,,(A,A',M',M,^,fT) (34) 

XfjX'fj.' cr j 

the factor 2 in Eq. 1341 accounts for the complex conjugate in Eq|2Sl since the matrix elements of 'L^ and L.S for 
parallel spins are imaginary and all are real. For a system we full translational symmetry and a single atom per 
unit cell this equation becomes 



<L;>^~2^ ^ ^(Aa|L;'|M'T)(MV|L.S|A"'a) (35) 

XfiX' fi' cr 

ief E'-E 



dE / dE' 

-oo J Ef 



X 

k 



Furthermore noting that [Liz"]xn = 2a{\a\L.S\^a) , Eq|^for < Lizn > can be transformed into: 



<U,„>^~4( J2 ^a(AT|L.S|A2t)(/i'T|L.S|A"'T)5]/y(A,A',Ai',M,'^,fT) (36) 

XfiX' fi' cr j 

For a system with full translational symmetry and one atom per unit cell, this yields 



< L,,„ >= -2^ J2 <^ t |L.S|a1 TXm' T |L.S|A"' T) (37) 



XjAX' jJ,' 



dE / dE' 

-oo J Ef 



X 

k 



^,Uxx' (k, S)X^v(k, S') + Mxx' (k, S)AA^,^(k, S') 



E'-E 



in which A/aa' (k, i?) = A/'AA'cr(k, i?). The generalization of this equation to systems with several atoms per unit 
cell is straightforward. 

It can be seen from Easl28l and 1361 that < Liz" > and the anisotropic part of AiJp-' are both given by quadratic 
functions of the direction cosines of the spin quantization axis relative to the crystal ones since the involved matrix 
elements of L.S are all proportional to one of these direction cosines (see Appendix A). These two functions present 

(2) 

some similarity, but spin- flip excitations contribute to AE^ but not to < Li^" >. However, if the exchange splitting 
is large enough compared to the d bandwidth, the spin up band is completely filled and the contribution of spin-flip 

(2) 

excitations to AE^ ' is negligible due to the large value of the energy denominator. In this condition, for each site i, 

(2) 

the anisotropy of AEl and < Li^'i > are proportional: 

AE^^\e,^) - Ai?f )(0,0) = -|(< Uz„{0,^) >-< L„,,(0,0) >) (38) 
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Note that this relation was aheady derived by Bruno P| for fee monolayers with a single atom per unit cell. Finally 
let us point out that for bulk cubic crystals with a single atom per unit cell both Ai?'^' and < Liz" > arc isotropic 
at the orders of perturbation considered. 

E. The ab-initio method 

For the sake of comparison we have also performed spin-polarized ab-initio calculations based on the Density 
Functional Theory (DFT) using the PWscf code of i^-ESFRESSO package |3| with ultrasoft pseudopotentials including 
non-linear core corrections. The calculations without spin-orbit coupling have been carried out within the GGA and the 
Perdew-Wang exchange-correlation parametrization, while the one including spin-orbit coupling have been performed 
within the LDA and the Perdew-Zunger exchange-correlation parametrization. The plane wave kinetic-energy cut-off 
was taken equal to 35Ry for the wavefunctions and 250Ry for the charge density and potential, which ensures a very 
good energy precision. 

F. Computational details 

When dealing with magnetic properties, and in particular magnetic anisotropy. the convergence of the total energy 
with respect to the number of fc-points has to be checked carefully. For all calculations involving magnetic anisotropy 
we checked that our results did not change by more than a few hundredth of meV (at most O.lmeV in the worst case). 
In the case of PWscf calculations the use of plane waves imposes a periodically repeated geometry and one must also 
avoid as much as possible electronic interactions by using large unit cells. The monatomic wires were separated by 
30a. u.. but for surfaces, we have been less demanding since we only calculated the magnetic moment and therefore 
the slabs were separated by approximately 17a. u. 

III. BULK MAGNETISM OF BCC AND FCC IRON 
A. Determination of the Stoner parameter 1^^/ from the magnetic transition in bcc iron 

In our TBHF model, the magnetism is entirely governed by the value of the Stoner parameter Idd'- It is well 
known that, in unsaturated magnetic materials like Fe, the magnetic moment is very sensitive to the precise value 
of the equilibrium interatomic distance. As a general trend, an expansion of the bulk lattice parameter leads to 
narrower (thus higher) density of states, which usually plays in favor of magnetism: it increases the magnetization 
in magnetic materials or it can trigger a magnetic transition in non-magnetic materials. Therefore a straightforward 
way to determine I^d' is to study the evolution of the magnetic moment as a function of the lattice parameter. In 
Fig. n the result of a series of TBLDA and TBGGA calculations on bulk bcc iron is shown for various values of the 
Stoner parameter. As expected the magnetic moment increases when the lattice is expanded but also when the Stoner 
parameter is increased. With Idd' = leV (TBLDA) and Idd' = l.lOey (TBGGA) we have been able to reproduce 
closely the results of PWscf calculations in a range of Wigner-Seitz radii (Rws) around equilibrium (the experimental 
bcc lattice parameter, 2.87A, corresponds to Rws ~ 2.67 a.u.). In the following we will keep these values fixed and 
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neglect any variation of these parameters with the local atomic environment. Let us however note that at very large 
lattice spacings the spin moment saturates (not seen in Fig^l at an "atomic" value of 2 /is for TBLDA and 4 /is for 
TBGGA. These two limits correspond to the different atomic configurations 3(PAs'^ and 3d^4s^ found in TBLDA and 
TBGGA, respectively, for a free Fe atom. Therefore the TBLDA gives a wrong atomic configuration which will have 
some consequences on the surface magnetism. Let us finally mention that in ab-initio calculations LDA and GGA 
yield very similar results as far as the magnetic moment is concerned. 




FIG. 1: Variation of the absolute value of the spin magnetic moment (per atom) of bcc Fe as a function of the Wigner-Seitz 
radius Rws for the values of the Stoner parameter I^^i (in eV) given on each curve. Left and right panels correspond to 
TBLDA and TBGGA calculations, respectively, compared to PWscf calculations in GGA. The dashed vertical line gives the 
experimental Wigner-Seitz radius at equilibrium. 



B. Fee iron 

The ground state phase of iron in normal temperature and pressure conditions is ferromagnetic (FM) bcc, at higher 
temperatures, the fee phase is stabilized but in a NM configuration. However it has been shown experimentally that 
thin films of iron can be stabilized in an fee structure |26l l27l |. This experimental work also showed the existence of 



various magnetic phases. It has also been known for a long time from theoretical works [23, l29l l30l l31l l32| that fee 
iron has a much more complicated magnetic structure than bcc iron. In the following, we present a study of the 
magnetic properties of fee Fe with TBGGA parameters. This will provide us with a first check of our model. 



1. Magnetic transition in fee iron. 

We have carried out a series of calculations on the fee phase of iron. Fig. |21 is similar to Fig. but for the FM 
(a) and antiferromagnetic (AFM)(b) bulk fee phases, the latter corresponding to a stacking of (001) planes in which 
spins of adjacent layers are opposite. In the FM case (Fig. EJi), the curves of appearance of a magnetic moment 
show a strong dependence on the magnetic moment Mq chosen as input to begin the self-consistency iterations (a 
similar behavior is also found with the PWscf code). For large Mq an abrupt transition from a NM configuration to 
a High Spin (HS) state occurs, for instance at Rws — 2.6a.u. when AIq ~ Sfis- For a small value of Mq a similar 
transition appears but at a much larger volume {Rws — 2.7a.u), while for intermediate Mq this transition is less steep. 
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The magnetic transition is much less abrupt in the antiferromagnetic phase (Figl^b) where two steps are observed: 
first, from a NM state to a low spin state (LS) and, second, from the LS state to the HS state. Note that a similar 
dependence on the input magnetic moment also exists in the AFM state (not shown on the graph). 




FIG. 2: Variation of the absolute value of the magnetic moment (per atom) for the ferromagnetic (FM) and antiferromagnetic 
(AFM) states of fee Fe as a function of the Wigner-Seitz radius Rws, for various input magnetic moments Mo in fiB obtained 
with TBGGA parameters. LS and HS denote respectively low and high spin states. 



2. Low and High spin ferromagnetic states: fixed spin moment calculations. 

The strong dependence on the input magnetic moment suggests the existence of metastable magnetic solutions. 
Therefore we have carried out TBGGA calculations using a fixed spin moment procedure for a series of Wigner-Seitz 
radii corresponding to the region of the magnetic transition. The behavior of the total energy as a function of the 
total moment M (Fig. |3J) reveals the existence of several local minima. In particular the curve at Rws = 2.67a.u. 
exhibits three minima (inset of Fig|3J): one at A'l — 0, one around M = 1.2/iB and one at M — 2.5/iB, corresponding 
to the NM, LS and HS states, respectively. This complex energy behavior is in agreement with the results obtained by 
Moruzzi et al. |3,'j | who showed for the first time the existence of three phases. It is clear that depending on the value 
of the input moment, the iteration loop will converge towards one of the three self-consistent (stable or metastable) 
magnetic states. 

C. Relative phase stability: comparison between TBLDA and TBGGA models. 

It is well known that DFT in the LDA predicts that, at low temperature, the fee AFM phase is the most stable one 
contrary to experiments [^Jl- However the right phase stability is recovered with the GGA. It is therefore interesting 
to investigate the ground state properties of bulk iron within our TB model. The results are presented in Fig. 01 it 
is found that TBLDA, similarly to DFT-LDA, gives the AFM state of fee iron as the most stable phase, while with 
TBGGA the FM bcc phase is found to be the ground state. These results are in perfect agreement with ab-initio 
findings, showing the ability of our model to reproduce rather complex magnetic behaviors. 
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FIG. 3: Variation of the TBGGA total energy per atom with the magnetic moment (fixed spin moment calculation for the 
ferromagnetic state) for fee Fe at several Wigner-Seitz radii Rws in a.u.. Note the presence of stable (or metastable) non 
magnetic (NM), Low spin (LS) and High spin (HS) states which is clearly seen in the inset. The zero of energy is arbitrary. 
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FIG. 4: Total energy per atom as a function of the Wigner Seitz radius Rws for ferromagnetic (FM), antiferromagnetic (AFM) 
and non magnetic (NM) states of bcc and fee iron. Panels a) and b) correspond to calculations performed with TBLDA and 
TBGGA, respectively. The zero of energy is arbitrary but is the same for all the curves in each panel. 

D. Influence of spin-orbit coupling 

The results presented above have been obtained without including spin-orbit coupling effects. The value of the spin- 
orbit coupling parameter ^ can be deduced from a comparison of the NM bcc band structure along a high symmetry 
direction of the Brillouin zone, for instance TH, obtained with our model and with the PWscf code. Indeed, the effect 
of spin-orbit coupling is to remove the degeneracy of degenerate levels when a matrix element of Hso exists between 
the corresponding eigenstates (see Appendix A). For instance it is easily seen that the six- fold degenerate level F^g/, 
corresponding to t2g spin-orbitals, are coupled by some matrix elements of Hso- Using the perturbation theory for 
degenerate levels, it is easily found that this level splits into a four- fold degenerate level at F25, — C/2 and a doubly 
degenerate level at Fjg/ -I- ^. From this splitting calculated with the PWscf code, we obtain ^ ~ 0.06ey. We have 
verified that with this value, our model is able to reproduce perfectly the spin-orbit coupling effects along TH. 

As seen in section 2, spin-orbit coupling is at the origin of the magneto-crystalline anisotropy. However, since this 
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anisotropy is of fourth order in ^, typical values for bulk materials are very small and of the order of 10~^ — 10~^ cV 
per atom for Fc. Co or Ni which makes the calculation of this anisotropy almost impossible, since it is beyond the 
accuracy of electronic structure methods. On the contrary, reliable values of the orbital moment, which is isotropic 
in the bulk to first order in perturbation, can be derived from Eg 1201 With our model and TBGGA parameters we 
find < >= O.Ol^iB in good agreement with experiments {O.OSfiB ISal)- 

As a conclusion, the TB results presented above (section 3.2 and 3.3) are in perfect agreement with DFT calculations, 
showing the ability of our model to reproduce rather complex magnetic behaviors. 



IV. (110) AND (001) SURFACES OF IRON 



We have then applied our TBGGA model to the study of the (001) and (110) surfaces of bcc iron. At the surface 
some atoms have a reduced coordination and therefore charge transfers as large as some tenths of electron are found 
in our model if the atomic levels £s, £p, Sd in Eq[21are kept at the values given by the MP equations. However, it 
is known that in metals, due to screening, the charge transfers are expected to be at least one order of magnitude 
smaller. To avoid unphysical charge transfers at surfaces the Hamiltonian is corrected by adding a term depending 
on the charge transfer SNi and of an average Coulomb integral U which must be large enough [U ~ 5eV) as shown 
in Eq|3 



A. Band structure of the (110) surface 



In a previous work on rhodium surfaces we showed that the charge quasi-neutrality is crucial to obtain a good 
description of the surface and resonant states [s^. Indeed these states are extremely sensitive to the energy shift 
induced by the rcnormalization of the intra-atomic terms of the TB Hamiltonian. Here we have carried out a TBGGA 
projected band structure calculation for the (110) surface of bcc Fe. The results are shown in Fig. [S] It can be seen 
that the position and size of pseudo-gaps in the band structure is significantly different for up and down spins. In 
particular along the TS and SH directions the pseudo-gaps are much larger in the minority spin band structure than 
in the majority spin one. As a consequence there are more minority spin than majority spin surface states. This is 
evidenced by the presence of a sharp down spin surface state around the Fermi level (indicated by an arrow in Fig. 
EJl which disappears in the up spin band structure. These results are in excellent agreement with previous ab-initio 
calculations [33 in particular, for the position and dispersion of the characteristic down spin surface state discussed 
above. 



B. Spin magnetic moments of Fe(llO) and Fe(OOl) surfaces 

It is well known that the lowering of coordination induces a narrowing of the density of states which usually enhances 
the magnetic moment. Consequently, it is expected that open surfaces should have larger surface magnetic moments 
than close-packed ones. We have therefore carried out self consistent TBGGA and PWscf GGA calculations for (001) 
and (110) surfaces. The (001) surface being more open than the (110) one, since each atom from their outermost 
layer looses 4 and 2 first nearest neighbors, respectively, we expect larger surface magnetic moments for the (001) 
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FIG. 5: TBGGA projected band structure for up (top) and down (middle) spins of a 20-layer (110) slab of bcc Fe with the 
experimental lattice parameter of 2.87A. The energy zero is the Fermi level. Surface or resonant states (i.e., states with more 
than 60% of their total weight on the first two outer layers) are represented in red and with thicker dots. A characteristic 
surface state of minority spin is indicated by an arrow. A schematic representation of the surface Brillouin zone and of the 
path in the reciprocal space is shown at the bottom. 

than for the (110) surface. Fig. shows that this general rule of thumb is well obeyed. Actually two clear features 
are seen in Fig. the magnetic moment is more reinforced on the (001) than on the (110) surface (+38% and +14%, 
respectively, for the outermost layer compared to the bulk ). Fricdel type oscillations arc present on the (001) surface 
while an almost monotonic decrease is obtained for the (110) surface. An excellent agreement is once again observed 
between PWscf and TBGGA results, in particular the spin moment is almost saturated on the outermost layer of the 
(001) surface in both calculations. Let us however note that the agreement is less perfect within TBLDA (not shown), 
which can be attributed to the wrong atomic configuration obtained in this model which deteriorates the spd charge 
distribution on the surface plane. 

C. Magneto-crystalline anisotropy 

For surfaces the magneto-crystalline anisotropy is usually one or two orders of magnitude larger than in the bulk. 
Indeed, it is well known that, contrary to the bulk, this anisotropy is of the second order in ^ at surfaces. Actually, 
we have seen in section 2.4 that second order perturbation theory predicts that the magneto-crystalline anisotropy is 
a quadratic function of the direction cosines (/ = sin^cosi/?, m = sm9smip,n = cos 9) of the spin quantization axis 
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FIG. 6: Variation of the spin magnetic moment (per atom) on successive atomic layers of (001) and (110) slabs (with 20 atomic 
layers) of bcc Fe obtained from the TBGGA model and PWscf code with GGA. Layer 1 corresponds to the outermost layer 
and layer 10 to a central layer. The value of the bulk magnetic moment is indicated as a reference. 

relative to crystal axes. By imposing the symmetry properties of the surface to this quadratic form the following laws 
are easily derived: 

AEf^ (6, ifi) ~ l^Ef^ (0, 0) = sin^ 9 (39) 

for the (001) surface with x and y crystal axes parallel to the edges of the square two dimensional cell, and: 

A£;f ^ (0, (^) - A^;f ) (0, 0) = i^l sin^ e + Af sin2 cos 2^ (40) 

on the (110) surface. For this surface the crystal axes are chosen as follows: the z axis is perpendicular to the surface 
and the y one is parallel to the second nearest neighbor direction in the surface. 

In Table^we present the results of TBGGA calculations on (001) and (110) slabs of two thicknesses. We have first 
checked the convergence of the magnetic anisotropy with respect to the number of fc-points in the first Brillouin zone. 
It is seen in Fig. [7|that a good convergence is obtained above 1000 fc-points. It is found that for both orientations the 
easy axis is perpendicular to the surface plane. The two surfaces, however, have a different behavior with respect to in- 
plane magnetization. For the (001) surface the energy is changing by only 4.10~^cV between a magnetic configuration 
along the square edge {lp = 0) and along the diagonal of the square [ip ~ 7r/4). This confirms our previous symmetry 
analysis which predicts no in-plane dependence of the energy at second order in perturbation. In the case of the (110) 
slab with 11 layers the in-plane energy variation is almost one order of magnitude larger and plays in favor of the 
second nearest neighbor atomic direction {ip — 7r/2). Consequently the anisotropy is very small in the yz plane. Let 
us also point out that isolated monolayers show a more pronounced out of plane anisotropy while the in-plane energy 
profile seems to be less corrugated than for thicker layers in the case of the (110) orientation. 
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FIG. 7; Convergence of the magnetic anisotropy -Etot(7i"/2, 0) — Etot{0,0) with respect to the number of k points Uk in the 
first BriUouin zone, for (110) and (001) Fe bcc, unsupported monolayers. The magnetic anisotropy is oscillating around its 
asymptotic value (full straight line) with an amplitude below ±0.02 meV when Uk is larger than 1000. 
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-0.17meV 



TABLE I: Magneto-crystalline anisotropy, Etot{0,ip) being the total energy (per surface unit cell) corresponding to a magneti- 
zation direction defined by the angles 9, if with respect to the crystal axes (see text) for slabs of bcc Fe with (001) and (110) 
orientations: n; is the number of layers and Uk is the number of k points in the first BriUouin zone used in the calculation. 



D. Orbital moment 



In presence of a surface the coordination is reduced and the symmetry is lowered leading to an enhancement of 
the spin moments as discussed above. At surfaces the orbital moment is also enhanced as demonstrated in several 
theoretical and experimental [ssl Is^ works. Our calculations on (001) and (110) slabs show that the component 

< Liz" > of the orbital moment on the magnetization direction is noticeably increased when i belongs to the outermost 
layer, especially for the (001) surface. On the second layer this component has almost recovered its bulk value for the 
(110) surface whereas oscillations occur for the (001) surface similarly to the behavior of the spin moment. Finally 

< Liz" >, contrary to the spin moment, depends sensitively on the direction of the magnetization with the same type 
of laws (Easl39l and l40|) as the magneto-crystalline anisotropy. However it is foimd that for 6 = tt/2 the ip dependence 
predicted by Eq2U]for the (110) surface is almost negligible. On the contrary the variation of < Liz" > with 6 is 
noticeable and, in this respect, the two surfaces behave differently: when the magnetization direction is rotated from 
6* = to 6* = 7r/2, < Liz" > decreases for the (001) surface while it increases for the (110) one. In addition, this 
variation is larger in absolute value on the (001) than on the (110) surface. 
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Finally it is expected that equation 1381 should be more obeyed for (001) surface than for the (110) since the spin 
moment at the outermost plane is saturated for the former and not for the latter. Assuming that the contribution to 
the magnetocrystalline energy comes from the outermost plane only, the ratio of the magnetocrystalline anisotropy 
to that of the orbital moment for surface atoms, has the wrong sign for (110), while for the (001) surface the sign 
is correct but the ratio is around three times smaller than ^/4. Indeed the exchange splitting of the (001) surface 
is smaller than the d bandwidth and the contribution of spin-flip excitations to the magnetocrystalline anisotropy 
cannot be neglected. 
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FIG. 8: Variation of the component of the orbital magnetic moment on the magnetization direction (per atom) as a function 
of the atomic layer in the (110) and (001) slabs (20 layers) for a magnetization perpendicular {9 — 0) or parallel to the surface 
(e = V2,^ = 0). 



V. STUDY OF THE MONATOMIC WIRE 

Although the study of the unsupported monatomic wire, i.e., a periodic linear chain of identical atoms with a 
single atom per unit cell, is somewhat academic it is however an interesting object for the following reasons: i) it 
can be used as a model since analytical TB results can be derived which arc useful for a theoretical understanding 
and a direct identification of the orbital character of the bands obtained in ab-initio calculations, ii) it also allows 
to investigate how a reduced dimensionality may modify magnetism in ferromagnetic metals or induce it in non- 
'QiE^i iii) last, but not least, such objects, several atom long, have been observed in break 
4J| (unfortunately not for Fe, Co or Ni). 

A. Non magnetic band structure of the monatomic wire 

The band structure of the non magnetic Fe monatomic wire, neglecting spin-orbit coupling, obtained from the 
TBLDA model and from the ab-initio PWscf code in the GGA are shown in FiglHlfor an interatomic distance of 4.29 
a.u.(2.27A), i.e., at the equilibrium distance predicted from the spin polarized PWscf GGA calculations. It can be 
seen that except for the upper band in the PWscf band structure, the agreement is satisfactory. 



magnetic materials |41 
junction experiments 
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Let us now identify the character of the bands by using the simplest TB model, i.e., in which overlap is neglected 
and hopping integrals are restricted to first nearest neighbors. Wc take the z axis along the chain. For a given spin 
the (9 X 9) hopping matrix can be rearranged into five square blocks on the diagonal: the first one involves s,^^ and 
d^^2_j.2 orbitals, the second and third ones are identical and involve {px, dzx) and {py, dyz) orbitals, respectively, finally 
the fourth and fifth ones are also identical and correspond to the d^y and d^^-y-^ orbitals, respectively. Consequently, 
taking into account spin degeneracy, the band structure consists in three two-fold degenerate bands of symmetry cr, 
two four-fold degenerate bands of symmetry tt and one four-fold degenerate band of symmetry 5. 
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FIG. 9: a) TBLDA and b) PWscf band structure of a non-magnetic Fe monatomic wire with an interatomic distance of 4.29 
a.u.. Each band is labelled by its symmetry character and degeneracy (including spin). 



The 5 band is easily identified as being the flattest band which disperses positively. Indeed its dispersion relation 
is that of a linear chain of d^y (or dj.2_j,2) orbitals with a small and negative hopping integral ddS. In Fe the p level 
is much higher in energy than the d level. As a consequence, the mixing between p and d orbitals is small and the 
TT bands split into a lower band with an almost pure dzx (or dyz ) character and a higher one with an almost pure px 
or py character. The first one shows a dispersion very similar to that of a linear chain of dzx (or dyz) orbitals, i.e., 
it disperses negatively since the corresponding integral ddi: is positive. The second one disperses positively since ppi: 
is negative, this band is present in the PWscf calculation, but is outside the energy range of Fig|51 in the TBLDA 
results. This means that the TB parameters relative to p bands are not very accurate but since the higher tt band is 
unoccupied, this inaccuracy will have no influence on our results for the ground state. The two remaining bands are 
the two lowest a bands. An analysis of the character of these bands using the TBLDA model shows that the lowest 
band has almost no p character, the weight of s and d^z^-r^ orbitals are almost the same at the F point while at the 
X point, the state is almost a pure d one. For the next a band the d^z^-r^ character decreases continuously from 
~ 0.5 to along TX, the pz character increases continuously from to 1 while the s character is ~ 0.5 at F, has a 
maximum at the midpoint and vanishes at the X point. 



B. Magnetic transition 




As in the bulk we have studied the appearance of a magnetic moment when the interatomic distance d increases 
using the TBLDA model as well as the PWscf code with GGA. As seen in Fig^ja low spin state (LS) is found at short 
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interatomic distances with a slowly increasing moment, then, around d = 3.8 a.u., the moment increases abruptly 
to reach a high spin state (HS) corresponding to the saturated configuration. The agreement between both methods 
is quite satisfactory, in particular, the HS state appears at the same interatomic distance. It should however be 
mentioned that the TBGGA fails to reproduce the wire band structure in the range of interatomic distances involved 
at the transition due to an incorrect position of the s level. Actually, for interatomic distances around 3.8 a.u. this s 
level is pushed at higher energies and instead of a LS/HS transition one finds a NM/HS transition. Let us recall that 
the TB parameters are fitted on bulk results only, and for interatomic distances and coordinence larger than in the 
wire. Therefore, the extrapolation of the law giving the atomic levels as a function of the atomic environment may 
fail for the wire. From the above study it appears that the TBLDA levels are much better than the GGA ones. As a 
consequence, in the following, LDA parameters will be used. 
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FIG. 10: TBLDA and PWscf spin moment of a monoatomic wire, as a function of the interatomic distance. 

C. Spin-orbit coupling effects 

1. Effects of spin- orbit coupling on the band structure from perturbation theory. 

The removal of degeneracies due to the spin-orbit coupling can be predicted by perturbation theory. As seen in 
FigEltwo types of degeneracies occur for the wire: in addition to spin degeneracy, either a band is degenerate for any 
value of k due to symmetry or the degeneracy is limited to band crossings. In the former case iJgo may produce a 
splitting of the bands while in the latter it may open a gap between the two crossing bands. 

Let us first consider the non magnetic case for which the band structure should be independent of 9 and ip. The 5 
band is four-fold degenerate. For 9 — 0(7r/2) Hso couples the dxy and d,j.2_y2 orbitals with the same (opposite) spins 
with a coupling matrix element (see Appendix A). Consequently the 5 band is unfolded symmetrically with a 
band splitting given by 2^. The same type of arguments applies to the lowest four-fold tt band if the small p character 
is neglected but the coupling matrix element is now ±i^/2 and therefore the band splitting is equal to ^. No removal of 
degeneracy is expected on the two- fold <t bands since there is no matrix elements of Hso between the d-i^^2_j.-2 orbitals 
of opposite spins. Finally in the unperturbed band structure (FigEJ there are two crossing points with the a and 5 
bands which are not coupled by Hso- Thus there is no removal of degeneracy. On the contrary the other crossing 
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points (between cr and n bands, or tt and 5 bands) are avoided. 

In the magnetic case the band structure is independent of ip for symmetry reasons but depends on the polar angle 
6 when spin-orbit coupling is taken into account. This is the consequence of the removal of spin degeneracy in the 
unperturbed state. Indeed for 6 — 0, Hgo couples states with the same spin and a splitting 2^(^) exists in the S 
(lowest tt) bands, for both spins, while for 6 = tt/2, the bands of up and down spins being well separated in energy, 
Hso has a negligible effect on the S and lowest tt bands. Finally, removals of degeneracy are expected at some crossing 
points. The number of opened gaps should be small at 6 — since, for bands of different symmetries, only states with 
opposite spins may be coupled by Hgo- As the exchange splitting is large the number of avoided crossings should be 
very small. On the contrary at 6* = 7r/2 a detailed analysis of the Hso matrix (see Appendix A) reveals that states 
with the same as well as different spins may be coupled and the number of avoided crossings is expected to be larger 
than for 9^0. 

The TBLDA calculation is presented in Fig. ^Jfor a magnetization parallel {9 = 0) or perpendicular {9 ~ 7r/2) to 
the wire. The results are in perfect agreement with the predictions of perturbation theory. Note that a rather good 
overall agreement is also obtained with PWscf calculations 45j. In particular the splitting of the 6 band in the latter 
calculations is 120meV, i.e., exactly 2^ which confirms that spin-orbit coupling is an intra-atomic effect and that ^ is 
a purely atomic quantity, i.e., it does not depend on the atomic environment. 
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FIG. 11: TBLDA band structure including spin-orbit coupling for a magnetic Fe monatomic wire (interatomic distance d=4.16 
a.u.) with a magnetization a) parallel and b) perpendicular to the wire. 



2. Magnetocrystalline anisotropy. 

Since the monatomic wire is a one dimensional system with the lowest possible coordination we expect a magne- 
tocrystalline anisotropy larger than at the surface. Moreover, due to the axial symmetry of the wire the energy will 
only depend on the angle 9 between the magnetization and the axis of the wire. By imposing the symmetry properties 
of the wire to the quadratic form in /, m, and n giving the magnetocrystalline anisotropy in second order perturbation 
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theory, we find: 

AEf ^ (6*, if) - AeI^^ (0, 0) ^ K[ sin^ 6 (41) 

In Fig. ^]we present the results of our TBLDA and PWscf calculations obtained at various interatomic distances. 
The variation of energy with 9 is perfectly fitted by the above equation. Interestingly an inversion of the easy axis 
is observed (i.e., a change of sign of K[). Indeed, for interatomic distances d larger than 3.78 a.u. (3.93 a.u. with 
PWscf) the easy axis is along the wire, while it is perpendicular when d is smaller. This inversion of the easy axis 
occurs at the LS/HS magnetic transition. A detailed analysis of the charge distribution shows that the magnetic 
transition is accompanied by a noticeable change of filling of the S bands. 
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FIG. 12: Variation of the magnetocrystalline anisotropy (per atom) of a monatomic wire as a function of the angle 6 from a 
TBLDA calculation (left panel) and PWscf calculation (right panel) for various interatomic distances d (in a.u.). 

3. Orbital moment. 

Let us first derive the expressions of the components of the orbital moment relative to the spin framework using 
perturbation theory. In a monatomic wire of atoms with only d orbitals, the band hamiltonian is reduced to five 
decoupled linear chains with a single orbital on each atom, thus JV\x'aO^, E) = JVx\aO^, E)S\x' and Eg 1361 becomes: 

< L; -2e^^(Aa|L;'|/ia)(Aa|L.S|/Ia)*/A^ (42) 

where /a^ denotes the term in Eg 1361 involving the summation over k but with A' = A and /i' = 

In order to find the angular dependence of < >, we use Ea ll8l with the matrix elements of and of Hso given 
in Appendix A taking = 0, for simplicity, since the wire has an axial symmetry. Then it is easily seen that: 

<L„. > = K,^..s\n2e (43) 
<L,." > = ^o."+^W'sin'^ (44) 
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while < ijj," >= as expected since y is perpendicular to the plane defined by the wire and the spin quantization 
axis. 

The results of the TBLDA calculations of < L,-^// > and < L^^u >, shown in Fig^jfor an interatomic distance 
d = 4.29 a.u., are in perfect agreement with this analysis. Note that < L^," > can reach a value as large as 
0.42/iB when the magnetization is along the wire and that the anisotropy of the orbital moment < L^^n {9 = n/2) > 
— < L^^ii{6 ~ 0) >= — 0.22/iB, is rather large. At the same interatomic distance (sec Fig ll2|) the corresponding 
magnetocrystalline anisotropy is 3.3.meV. Thus the ratio of this anisotropy to that of the orbital moment is equal to 
-15meV which is in perfect agreement with Eg 1381 f— f /4 ~ —ISimeV). This was rather expected since the spin up d 
bands are completely filled and well separated (~ 3eV) from the spin down ones at this interatomic distance, both 
bands being much narrower than at the surface. 
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FIG. 13: Variation of the components of the orbital moment relative to the spin framework for a monatomic wire (interatomic 
distance d — 4.29 a.u.) as a function of the direction of magnetization given by the angle 6 obtained with TBLDA parameters. 

At shorter distances for which the LS state is found, Eq|3Hlis no longer valid and even does not predict the sign of 
this ratio since the magnetocrystalline anisotropy changes sign, contrary to that of the orbital moment. 

Finally let us note that in the presence of spin-orbit coupling the spin moment is very slightly anisotropic: typically 
< Sz" > varies by about 0.02/iB between 9 — and 6 = tt/2. Moreover, similarly to the orbital moment, a small 
component < Sx^' > is present when e]0,7r/2[. 



VI. CONCLUSION 



We have shown that starting from the parametrized spd TB model set up by Mehl and Papaconstantopoulos 
|2l| . adding a Stoncr-like spin polarization term and a spin-orbit coupling term, i.c, introducing only two additional 
parameters 7^^' and ^, we have been able to describe in detail the magnetic properties (spin moments, orbital moments 
and magnetocrystalline anisotropy energy) of iron in systems of various dimensionalities and coordinences. Whenever 
possible the results have been compared with those of the PWscf code or other existing theoretical or experimental 
data, and the agreement is excellent. Our simple TB model allowed us to derive some new results. For example 
in the case of surfaces we have studied the variation of the component of the orbital magnetic moment on the spin 
quantization axis < > as a function of depth, and shown that the enhancement of this quantity is still noticeable 
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on the second layer for the (001) surface but almost cancels on the third and innermost layers. In the wire wc have 
found that the easy axis of magnetization is along the wire at the theoretical equilibrium distance but can switch to 
the perpendicular direction under compression. Moreover < L^" > is strongly enhanced and highly anisotropic since 
for z" along the wire it is about 0.42/is and decreases to half this value for the perpendicular direction. In addition 
for intermediate orientations the orbital moment has a non-negligible component < L-c" > perpendicular to z" and 
in the plane made by the wire and z" . This component is also anisotropic and proportional to sin 26 as predicted 
by perturbation theory. Finally we have shown that the law of Bruno Q relating the MAE and the anisotropy of the 
orbital moment, is perfectly obeyed for the wire at equilibrium distance, which is not the case for the surfaces, even 
for the (001) orientation at which the magnetic moment is saturated. 

The success of our model opens up the possibility of obtaining accurate results on other elements and systems with 
much more complex geometries. 
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APPENDIX A: THE MATRIX ELEMENTS OF THE ORBITAL MOMENT AND SPIN-ORBIT 

COUPLING OPERATORS 

The matrix elements (Acr|L.S|/i(T'), where A and /I are the angular parts of d orbitals {dxy, dy^, dzx, d^^-y^, d^z^_r2) 
centered at the same site as the angular momentum operator, can be easily calculated as a function of the direction 
(defined by the polar and azimuthal angles 9 and ip relative to the crystal axes) of the spin quantization axis z" which 
is often taken in the magnetization direction that may be different from the z axis of the crystal. This is achieved by 
choosing a new coordinate framework x'", y", z" referred to as the spin framework which is obtained by rotating first 
the xy axes of the crystal by the angle ip around z. This gives a new framework x', y', z' which is then rotated by the 
angle around y' . Using the matrix elements of L in the basis of real d orbitals: 
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one finds: 



(Act|Lz|/1(t) = 



/ 2i o\ 

i 



-2i 

\ 0/ 



(A T IL.SI/2 T) = 



/ 

^ sin sin ( 
I cos </? sin 6 
— i cos 6* 




2 sm sm t 




-J- cos^ 



5 COS If sm t/ 
cos Lp sin ( 



cos ip sm y I cos y 

I cos ^ cos sin ( 

^ sin iy9 sin ^ 

I sin sin 9 

sin 1^ sin 



— cos (/3 sin t 
^ sin </? sin 9 





(AT IL.SI/ii) 



/ 

i sin 9 

V 



—i sin 61 



where f{cp,9) = cos(^ + ismipcos9 and g{ip,9) ~ s'mip — icos(y9Cos6'. The other blocks of the (10 x 10) spin-orbit 
matrix are obtained from the relations: 



(Aj |L.S|/2t) = -(AT |L.S|Mi)* (Al) 
(AT IL.SImT) - (AT |L.S|a2T)* (A2) 

in which * denotes the complex conjugate. In addition a very useful relation has been derived by Bruno 



Re[(A T |L.S|m i){fl' I |L.S|A' T)] + (A T |L.S|/1 TXm T |L.S|A' T) = Cst. 



(A3) 



Let us note however that the spin quantization axis could have been taken along the z axis of the crystal, in which 
case the spin-orbit matrix elements would be given by the above matrices with 9 = ip = but the spin polarized term 
of the Hamiltonian becomes: 



cos9 exp{—i(p) sm9 
exp(i(^) sin 9 — cos 9 



where (Xi means the direct product of matrices. 
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In the case of collincar spins and in the presence of spin-orbit coupling the first point of view is more convenient 
to treat spin-orbit couphng effects within perturbation theory since _ffso describes the perturbation completely. The 
second point of view is preferable when dealing with non-coUinear spins. Indeed this avoids the transformation of the 
inter-atomic part of Htbhf since, in that case, the spin fimctions arc not the same at the two sites. 

APPENDIX B: CALCULATION OF THE ORBITAL MOMENT FOR A NON-ORTHOGONAL BASIS SET 

Let us generalize EqElto take overlap into account. If we note that the integral in Eg 1131 gives the population of 
the spin-orbital \ilm.a), an obvious generalization is to replace this population by the MuUiken one. Thus pumaiE) 
becomes: 

pum.{E) = Re 5] arz™.«?;'™v4™' - i^n) (Bl) 

i I 7n 
n 

with: 



\i^n) - al^^\ilma > (B2) 

ilrntr 

and 



girn,i'm' ^ l^Uma\i'l'm'a) (B3) 



thus 



<L,,. >=Re ^<L.<n'.,A^'''^' (B4) 

n occ 

and, after simple algebraic manipulations: 



Im.i I m a 
h OCC 

The generalization of Eg 1161 yields: 



<L;'>=Re Y (V■nNW)[L;']^„,^,^„,45-l]'^;';"''''''"'(^'^'m'a|V'„) (B6) 

Im^l" m" ^i'V m' (7 
n occ 

and, in the basis IzAcr): 



<L:'>=Re Y (V'«NAa)[L;']A^[5-i],i:r(zVa|V'„) (B7) 

n OCC 
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with: 



(B8) 



so that we get finally: 



<L:' >=Rc ^ <:tx[L;']AM5^rc?,, (B9) 



Xfii' ua 
n occ 
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